Estudiante Ing. Forestal
NASA
“Desde 1880, el nivel del mar global ha aumentado 20 centimetros. Para el 2100 se proyecta que aumente entre 30 y 122 centimetros más.”
Fuentes
pacman::p_load(TSstudio, tsbox, tsoutliers, tidyverse, equatiomatic,
lubridate, xts, magrittr, forecast, dygraphs)
url <- 'https://gml.noaa.gov/webdata/ccgg/trends/co2/co2_mm_mlo.txt'
cov <- read.table(file = url) %>%
select(3, 4) %>%
rename(fecha = 1, co2 = 2) %>%
mutate(fecha = date_decimal(fecha) %>%
floor_date(unit = 'month')) %>%
filter(fecha <= '2013-12-01')
data <- read.csv('gmsl.csv') %>%
mutate(fecha = date_decimal(Time) %>%
floor_date(unit = 'month')) %>%
rename(gsml = 2) %>%
select(4, 2)%>%
filter(fecha > '1958-02-01') %>%
left_join(y = cov) %>%
{xts(x = .[,2:3], order.by = .$fecha)}
rm(cov)
dygraph(data, main = 'Nivel medio del mar global vs
Concentración de CO2 en la atmósfera', height = 400) %>%
dySeries("co2", axis = 'y2', label = 'CO2') %>%
dySeries('gsml', label = 'GSML') %>%
dyAxis('y', label = 'Nivel del mar (mm)') %>%
dyAxis('y2', label = 'CO2 (ppm)', independentTicks = T) %>%
dyRangeSelector()
z test of coefficients:
Estimate Std. Error z value Pr(>|z|)
ar1 0.797070 0.054991 14.4947 < 2.2e-16 ***
ar2 -0.095533 0.054919 -1.7395 0.08194 .
ma1 -0.339090 0.040237 -8.4274 < 2.2e-16 ***
ma2 0.208116 0.038717 5.3753 7.646e-08 ***
ma3 -0.693544 0.028246 -24.5535 < 2.2e-16 ***
drift 0.208518 0.041471 5.0281 4.955e-07 ***
co2 -0.171228 0.087716 -1.9521 0.05093 .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AIC: 2595.726
z test of coefficients:
Estimate Std. Error z value Pr(>|z|)
ar1 0.794860 0.055187 14.4031 < 2.2e-16 ***
ar2 -0.092967 0.055205 -1.6840 0.09218 .
ma1 -0.341537 0.040571 -8.4182 < 2.2e-16 ***
ma2 0.204242 0.039030 5.2329 1.668e-07 ***
ma3 -0.692204 0.028275 -24.4809 < 2.2e-16 ***
drift 0.187589 0.039109 4.7965 1.614e-06 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AIC: 2597.526
Holt's method
Call:
holt(y = train[, "gsml"], h = 12, initial = "optimal")
Smoothing parameters:
alpha = 0.9999
beta = 0.6589
Initial states:
l = -49.0075
b = -2.4967
sigma: 2.4035
AIC AICc BIC
5429.938 5430.030 5452.384
\[ \begin{alignat}{2} &\text{let}\quad &&y_{t} = \operatorname{y}_{\operatorname{0}} +0.21\operatorname{t} -0.17\operatorname{co2}_{\operatorname{t}} +\eta_{t} \\ &\text{where}\quad &&(1 -0.8\operatorname{B} +0.1\operatorname{B}^{\operatorname{2}} )\ (1 - \operatorname{B}) \eta_{t} \\ & &&= (1 -0.34\operatorname{B} +0.21\operatorname{B}^{\operatorname{2}} -0.69\operatorname{B}^{\operatorname{3}} )\ \varepsilon_{t} \\ &\text{where}\quad &&\varepsilon_{t} \sim{WN(0, \sigma^{2})} \end{alignat} \]
\[( 1 \ - \ 0.79B \ + \ 0.09B^2) \ \cdot \ (1 - B) \ \cdot \ (y_t - 0.19t) \\ = (1 \ - \ 0.34B \ + \ 0.2B^2 \ - \ 0.69B^3) \ \cdot \ \epsilon_t\]
Series:
Regression with ARIMA(2,1,3) errors
Coefficients:
ar1 ar2 ma1 ma2 ma3 drift xreg
0.7971 -0.0955 -0.3391 0.2081 -0.6935 0.2085 -0.1712
s.e. 0.0550 0.0549 0.0402 0.0387 0.0282 0.0415 0.0877
sigma^2 = 2.993: log likelihood = -1289.86
AIC=2595.73 AICc=2595.95 BIC=2631.63
No outliers were detected.
prediccion <- forecast(object = modelo, h = 12, level = 95, xreg = test[,'co2'])
serie <- prediccion %>%
as_tibble() %>%
ts(start = c(2013, 1), end = c(2013, 12), frequency = 12) %>%
{cbind(data[, 'gsml'], .)}
colnames(serie) <- c('gsml', 'predicho', 'lwr', 'upr')
serie %>%
dygraph(height = 350) %>%
dySeries('gsml', label = 'Observada') %>%
dySeries(c('lwr', 'predicho', 'upr'), label = 'Predicho') %>%
dyRangeSelector()serie <- holt %>%
as_tibble() %>%
select(1, 4, 5) %>%
rename(predico = 1, lwr = 2, upr = 3) %>%
ts(start = c(2013, 1), end = c(2013, 12), frequency = 12) %>%
{cbind(data[, 'gsml'], .)}
colnames(serie) <- c('gsml', 'predicho', 'lwr', 'upr')
serie %>%
dygraph(height = 350) %>%
dySeries('gsml', label = 'Observada') %>%
dySeries(c('lwr', 'predicho', 'upr'), label = 'Predicho') %>%
dyRangeSelector()Series de Tiempo Univaridas